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1 Introduction 



In current semiconductor technology, vias and interconnect leads are manufac- 
tured by physical vapour deposition (PVD, sputtering) of a copper seed layer 
and subsequent electrochemical metal deposition. Steadily decreasing feature 
sizes reduce the lead's diameter and thus increase the resistivity Furthermore, 
current density increases, which promotes electromigration of copper ions. Also, 



vias with high aspect ratios are hard to fill by sputtering. Figures 1.1 and 1.2 
show an schematic representation of such structures. 

In the strive for a sustainable metallization technology. Atomic Layer Deposi- 
tion (ALD) has recently been proven as a means to grow copper oxide layers on 
spatially structured as well as on large plane substrates [1|. Thomas Wachtler 
presented a hybrid approach to divide the metallization approach into an cop- 
per oxide ALD and a subsequent reduction step. Direct deposition of metallic 
copper was not possible in the desired quality because copper tends to form 
islands and clusters (agglomeration), but the field is also actively investigated 

|2il- 

ALD is a promising technology because it is known to produce void free 
and homogeneous films. Benefits are expected since it is possible to deposit a 
seed layer in high-aspect-ratio vias by ALD and because films can be deposited 
more homogeneously than with current technology. A film with less voids will 
exhibit lower resistivity and electromigration. The major obstacle in employing 
this technology is the reduction of copper oxide to metallic copper, which has 
to be carried out at low temperatures to avoid agglomeration of the resulting 
copper film. 




Figure 1.1: SEM cross-section of AMD Opteron™ and AMD Athlon™64 mi- 
croprocessors, showing the 9-metal interconnects hierarchy (lEDM 
2003, Nanofair 2003), Photo courtesy AMD Saxony (now GLOBAL- 
FOUNDRIES) - Dresden, Germany 
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Figure 1.2: Schematic illustration of the metallization layers in recent semicon- 
ductor devices, picture courtesy of Steve Miiller. 



The actual reduction reaction on the samples is hard to characterize experi- 
mentally because the grown film is not epitaxial and contains both copper(I)- 
and copper(II)-oxide. Moreover, measurements would have to be performed 
without exposing the samples to ambient conditions after the ALD process, and 
spectroscopic equipment is currently not available at the ALD chamber. In cur- 
rent investigations samples were exposed to air during transport between the 
ALD chamber and the spectroscopic measurements. Studying the most basic 
reactions in a theoretical work at the electronic structure level is thus a first step 
towards a more thorough understanding of the process. 

In order to gain a deeper understanding of the reactions involved from a 
theoretical point of view it is necessary to find a valid model system for the 
copper oxide surface that can be used to investigate reaction mechanisms and 
test proposed reactions for their energetic profile. Hence, the first task in this 
work will be to select a feasible copper oxide structure for electronic structure 
calculations. Due to the accuracy requirements and the system size of about 100 
atoms, density functional theory (DFT) appears as the optimal choice for all cal- 
culations as it delivers a reasonable balance between computational complexity 



and accuracy. This is discussed in more detail in section [Z2 



Formic acid and its active decomposition products, hydrogen and carbon 
monoxide, are known to act as reduction agents for copper oxide. Decomposi- 
tion of formic acid and reactions with different adsorbed species on the surface 
can be investigated once a surface model is available. 



1.1 Preliminary Work 

During a previous internship at Fraunhofer ENAS the unimolecular thermal 
decomposition of formic acid in gas phase was modelled |j3|]. Two possible 
decomposition paths with similar reaction barriers and enthalpies were found: 

HCOOH -^ CO2 + H2 A£ = -3.5 r (1.1) 

mol 

HCOOH ^C0 + H2O AE = 6.3 ^ (1.2) 

mol 

Reaction barriers were estimated to 70 kcal mol"^ for both reactions. Reaction 
energies above are from [4|. Essentially the conclusion was, that under the 
reactor conditions described in Thomas Wachtlers work [IJ, formic acid will not 
dissociate fast enough to supply noticeable amounts of CO or H2 to the sample. 
The introduction to chemical kinetics in the Eyring model is a result of that work 
(See section [Z6l >. 



1 .2 Known Reactions and Issues 

Several experiments on reduction of copper oxide have been performed by 
Miiller et al. [5]. It has been shown that reduction of a CuO film to metallic 
copper can be done by hydrogen radicals on a ruthenium or nickel substrate, 
a reduction to CuiO is possible by molecular hydrogen on these substrates. 
On a cobalt substrate reduction proceeds to CU2O with hydrogen radicals. All 
reductions have been carried out below 300°C to avoid agglomeration of the 
deposited film. Complete Reduction by molecular hydrogen is possible at higher 
temperatures, but then agglomeration becomes an issue. 

Using formic acid as a reduction agent, Wachtler and Miiller observed no 
reaction on a tantalum nitride substrate, but complete reduction to metallic cop- 
per took place on a ruthenium substrate [5|. Suspecting a catalytic influence of 
ruthenium, films were prepared with a mixed copper and ruthenium precursor 
on TaN. It was shown, that the films contained 1 to 5 percent of ruthenium and 
that all CuO fractions could be reduced by formic acid. To date, it is not clear 
whether remaining CU2O fractions are due to incomplete reduction or due to 
re-oxidation during transport of the samples. 

From other experiments it is known, that clean Ruthenium and Nickel surfaces 
catalytically promote the decomposition of formic acid ||6|. 

Poulston et al. performed temperature-programmed desorption studies of 
formic acid decomposition on both CuO and CU2O [7]. They report a formate 
adsorption of formic acid, formate decomposition to gaseous CO2 and adsorbed 
hydrogen as well as desorption of molecular hydrogen and water formed from 
adsorbed hydrogen and lattice oxygen. Formic acid was adsorbed at 300 K and 
two desorption peaks of all species occur at 430 K and 545 K, where the latter 
one is more intense. 

Reduction of CU2O with CO and H or H2 is already well understood at high 
temperatures because of its application in automotive exhaust catalysts. 

The reduction of copper oxides with hydrogen and carbon monoxide has 
been investigated by Kim, Wang, et al. [Sl |9i. In a temperature programmed 
reduction using carbon monoxide as reduction agent, the reduction of CuO was 
found to proceed either directly to metallic copper when high amounts of CO 
were supplied or via formation of CU2O when CO supply was limited. In both 
cases a temperature dependant induction period was observed. The reduction 
of CU2O was reported to proceed slower than the reduction of CuO. Similarly, 
in a constant-temperature reduction using molecular hydrogen as a reduction 
agent, an induction period involving the embedding of hydrogen into the bulk 
oxide was reported from in-situ time-resolved x-ray diffraction. The activation 
barrier was, again, higher for CU2O (27.4 kcal mol"^) than for CuO (14.5 kcal 
mol"^). Additionally, density functional calculations were carried out on a bulk 



model to explain the transition from CuO to CU2O instead of CU4O3. 

Goldstein and Mitchell have recently measured reaction rates of copper oxide 
reduction with carbon monoxide using copper(I) and (II) oxide powder in a 
pressurized thermogravimetric analyser pO| . They report an activation energy 
of 20 kj mol"^ and 25 kj mol"^ for CuO and CU2O, respectively 



On Ruthenium, a reaction model was proposed by Sun and Weinberg |11| 
where two formate intermediates stabilize on neighbouring sites and a "hot 
hydrogen" breaks the C-O bond of a second formate, resulting in one desorbed 
CO2, one adsorbed hydroxyl(-OH) and one adsorbed formyl(-CHO). The ad- 
sorbed species react to form adsorbed CO and either desorbed H2O or desorbed 



H2 and adsorbed oxygen. Figure 1.3 illustrates the mechanism 
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Figure 1.3: Reaction model for bimolecular decomposition of formate on Ru, as 
proposed by Sun and Weinberg |[lT|. Figure reprinted from |6l. 



Copper oxide surface structures have been investigated by several groups 
with a variety of experimental and theoretical methods. 

Schulz and Cox performed photoemission and low-energy-electron-diffraction 



(LEED) studies of CU2O (111) and (100) surfaces |12|. The (111) surface was re- 
ported to stabilize in a stoichiometric (1x1) form after ion bombardment and 



annealing in vacuum. A ( V3 X V3)R30° periodicity was also observed during 
studies of catalytic chemistry due to a ^ monolayer of oxygen vacancies. The 
(100) surface exhibited four different periodicities in the LEED studies. A recon- 
structed, copper terminated surface with a (3 V2 x V2)K45° periodicity could 
be prepared by ion bombardment and vacuum annealing. After high oxygen 
exposures an unreconstructed, (1 x 1), oxygen terminated surface was reported. 
A DFT study by Soon et al. predicted a Gibbs free energy preference for a 
Cu-lean CU2O (111) surface where the Cucus species are vacant above 300 K and 
under oxygen exposure p3|. Figure 1.4 shows the calculated free energies (a). 



an optimized bulk terminated surface structure (b and c) and the low energy 
structure with Cucus vacancies (d). 

Another experimental study by Onsten et al. using LEED and scanning tun- 



nelling microscopy (STM) resulted in two models for the (111) surface |14|. In 
the first model the (1 x 1) form is the ideal bulk terminated surface and the 
( V3 X V3)R30° form is due to oxygen vacancies, in agreement with the re- 
sults from Schulz and Cox. The second model agrees with the first one, but 
Cucus species are vacant in both surface forms, which is supported by the find- 
ings of Soon et al. 

The polar CU2O (100) Cu"^ terminated surface was investigated by Nygren et 



al. in a theoretical model |15|. The (1 x 1) reconstruction with half of the Cu"^ 
placed at the opposite site of the crystal gave the lowest surface energy. The 
(3 V2 X V2)R45° reconstruction was not stable within their model. 
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Figure 1.4: A plots the surface free energy for the various Cu20(lll) surface 
terminations as a function of oxygen chemical potential. B and C 
show the top and side view of the optimized atomic structure of 
Cu20(lll), respectively. D shows the top view of the low-energy 
structure without Cucus species. 

The upper Cu atoms are shown as large orange spheres, with the 
oxygen atoms denoted by small red spheres. The labels "CSA" and 
"CUS" stand for coordinatively saturated, and coordinatively unsat- 
urated, respectively. Note that Ocus arid Ocsa will be called Osuf and 
OsuB in this work, respectively. Figure reprinted from llSjl. 



1 .3 Overview of Reactions and Species involved in 
Formic Acid Decomposition 

In the following section a list of reactions is presented that an formic acid 
molecule might undergo. This list is not exhaustive, other reactions may oc- 
cur, but the most likely reactions are listed and may be used as a starting point 
for later investigations once a copper oxide surface model is available. 

Formic acid may decompose thermally prior to adsorption : 

HCOOH -^ CO2 + H2 (1.3) 

HCOOH -^CO + H2O (1.4) 

or adsorb and then proceed to possible decomposition products: 

HCOOH(g) -^ HCOOH(,,3) (1.5) 

HCOOH(,,3) -^ H(,,^) + HCOO(,,3) (1.6) 

-^ OH(,,,) + HCO(,,3) (1.7) 

-^ CO(,,3) + H20(,,3, (1.8) 

The reaction mechanism for formic acid decomposition on Ru proposed by Sun 
and Weinberg may as well be existent on CuxO 

2 (HCOO(,,3)) -^ C02(g) + OH(,,3) + HCO(,,3) (1.9) 

OH(,,3) + HCO(,,3) -^ CO(,,3) + H20(g) (1.10) 

-^ CO(,,3) + 0(^,3) + H2(g, (1.11) 

Different reactions with previously adsorbed species are possible 

HCOOH(,,3) + H(^,3) -^ HCOO(,,3) + H,^^^ (1.12) 

-^ HCO(,,3) + H20(,,3, (1.13) 

HCOOH(,,3) + OH(,,3) -^ HCOO(,,3) + H^O^,,^^ (1.14) 

-^ COOH(,,3) + H20(,,3^ (1.15) 

as well as reactions with surface oxygen and adsorbed species 

OsuF-HOOCH + H(^,3) -^ HCOO(,,3) + H,0^^^ (1.16) 

A few of these reactions will be modelled to test the model system for its 
applicability. 



2 Theoretical Background 



Electronic structures are evaluated by means of density functional theory, which 
is considered the standard method for simulations at this scale because it allows 
easy (in terms of computer time) treatment of systems of up to some hundreds 
of atoms at reasonable accuracy Other, wave function based, techniques are 
available, where MP2 has comparable accuracy and is applicable to systems 
of 50 to 100 atoms. High accuracy methods like CCSD(T) and CI reach their 
limits at about 20 atoms. In order to predict possible reaction paths and to 
estimate reaction rates, reactants and products of a reaction are modelled at 
atomic scale and structures are optimized. Structure optimization means to find 
the atomic coordinates that exhibit a locally minimal DFT energy Eq. Comparing 
the reactants' and products' energy sums yields the reaction energy AE. 

AE= Yj Eo,p- Yj ^0'^ (^-^^ 

products reactants 

Vibrational effects are accounted for in a harmonic oscillator model, yielding 
respective enthalpies. Hq is the enthalpy that includes only the zero point 
vibrational energy (ZPVE), when higher vibrational states become occupied the 
enthalpy will be a function of temperature H{T). 

The reaction enthalpy will only provide insight as to whether the reaction is 
endothermic or exothermic. Estimating reaction rates requires knowledge of 
the activation energy. It is defined by the structure (and thus relative energy ) 
of the reaction's intermediate structure (transition state, TS). A transition state 
is characterized by a saddle point on the energy hypersurface. 

Throughout this work, AE and AH will be reaction energies and enthalpies, 
Ea will be an activation energy. Ha denotes an activation enthalpy evaluated 
between the zero point vibrational energies and Ha(T) the activation enthalpy at 



finite temperature. See Figure 2.1 for an illustration of the various energy and 
enthalpy terms. 

The following chapter shall first introduce the path from a quantum mechan- 
ical many body problem to the basic DFT formalism. Then, this technique is 
applied to geometry optimizations, transition state searches and the transition 
state model for chemical kinetics. 
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Figure 2.1: Energy profile of a hypothetical reaction, containing definitions of the 
various energy and enthalpy terms. Parabolas illustrate an harmonic 
approximation of the potential around the structure with discrete 
vibrational states. 
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2.1 The Schrodinger-Equation 

In quantum theory, the fundamental problem is to find the many-body wave 
function \p of the system. One can then evaluate expectation values of physical 
properties by applying the corresponding operator to the wave function. To 
find the wave function, the energy eigenvalue problem (ScHRODiNGER-Equation) 
must be solved. ScHRODiNGER's-Equation in time-independent form is: 

H|i/') = E\4^) (2.2) 

where ip is the wave function and H is the HAMiLxoN-operator. The Hamilton 
Operator is the total energy operator and for this problem will be the sum of 
a kinetic energy operator T of electrons and nuclei and potential energy opera- 
tors V that describe the electron-nucleon interaction and the mutual interaction 
among either electrons or nuclei. In atomic units (a.u.) these are: 

^ = f+V (2.3) 

H = te + t, + yne + t/ee + t>nn (2.4) 

N A4 N M ry 

= _ y lv2 - y J_vi - y y ^ 

L^2 ' L-L 2Ma ^ ^^nA 

i=l A=l "^ i=l A=l ^^ 

N N . M M „ ry 

!=1 i>i ') A=l B>A ^^ 

For electrons and nuclei indices {i,A) and position vectors {f[,Rp) are in lower 
and upper case, respectively. The first two terms cover the kinetic energy of 
electrons and nuclei, where V^ and V^ are the Laplacian differentiation oper- 
ators with respect to the coordinates of the zth electron and Ath nucleus. Za 
is the nuclear charge of atom A. The latter three terms represent the coulomb 
interaction between electrons and nuclei, electron pairs and nucleus pairs where 
'^■iy ''iA/ ''ab are the mutual distances between two electrons, one electron and one 
nucleus and between two nuclei. 



One common simplification is the Born-Oppenheimer- Approximation |16|: 

the set of spatial variables is divided into electronic coordinates r and nuclear 

-» 
coordinates R and the wave function is expressed as the product of two separate 

ones for electrons and nuclei 

i/^(f) = i/;e(r,R)-i/^n(^) (2.6) 

Applying this separation to the Schrodinger-Equation results in equations for 



11 



electrons 

H-4>e{r,R) = E,-,p,{f',R) (2.7) 

(te + tn + 1/ne + t^ee) " M^, R) ^ E, ■ l/^eC^ R) (2.8) 

and nuclei 

H-iPn{R) = E,-iPn{R) (2.9) 

(tn + V^) ■ XPniR) = En " M^) (2-10) 

which are coupled via the T^ ■ \pg{r,R) term. Decoupling the equations by 
eliminating this term is known as the Bom-Oppenheimer approximation. The 
electronic wave function will then depend parametrically on the nuclear coor- 
dinates ( ie. i/^e = 4'ei^' [^]) )' ^^t they will not count as degrees of freedom in the 
system. This introduces the concept of a potential energy surface on which local 



minima represent stable structures (see Section 2.5 ). 



The BoRN-OppENHEiMER-Approximation usually works very good. Jensen 



1 17 1 argues that it introduces an error of the order of 10 au which is usually 



negligible when compared to other errors. 

A more thorough and rigorous discussion can be found (among others) in 1 17 



18 19|. 



2.2 Density Functional Theory 

Contrary to the wave function approach, density functional theory focuses on 
the electron density p instead of the wave function: 

p(r-) = |^(r-)p = ^(r-)-i/.*('1 (2.11) 

The total energy E can then be expressed as a functional of electron density, 
where ENe is the energy due to nucleus-electron interaction, T is the kinetic 
energy of the electrons, / is the electronic coulomb interaction and E^a covers all 
non-classical effects (self-interaction correction, exchange, correlation). 



E = E[p] 
= T[p] + J[p] + E^,[p] + Enci[p] (2.12) 

Hohenberg and Kohn pO| have proven that the ground state electronic energy 
is given by a unique functional of the true ground state electron density and that 



12 



the ground state electron density minimizes the total energy of the system. The 
electron density is a function only of the spatial variables, instead of all electron 
coordinates (spatial and spin), which drastically reduces the complexity of the 



problem. However, modem DFT codes, by suggestion of Kohn and Sham |21 1, 
reintroduce these degrees of freedom by generating the density from a basis set 
of atomic orbitals (linear combination of atomic orbitals, LCAO). 

In Kohn-Sham Theory an auxiliary non-interacting reference system is con- 
structed: 



- N N 



(2.13) 



where Vg is an effective local potential. This has the advantage that the 
main contribution to the real, interacting systems energy (kinetic energy) can be 
calculated exactly. The corresponding ground state wave function can then be 
expressed as a Slater determinant: 



0s = 



/N! 



<Pi(4) (piixi) ■■■ (pNixi) 

(PliXi) (pi{X2) ■■■ (Pn(4) 

<Pi(^n) <P2(^n) ••• <Pn(xn) 



(2.14) 



where the cpi are the eigenstates of the one electron Kohn-Sham operator f^^ and 
Vc, is chosen such that the resulting density equals that of the real, interacting 
system: 



N 



(2.15) 
(2.16) 

(2.17) 



The energy of a non-interacting system can now be calculated exactly, but only 
the energy of the interacting system is physically meaningful. The (interacting) 
electronic energy functional is divided as follows: 

f[p] = rs[p] + /[p] + Exc[p] (2.18) 

where Exc is the exchange-correlation energy, which, by comparison with equa- 



tion 2.12 is defined as 



Exc[p] = T[p\ - TAp\ + End[p]. 



(2.19) 
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Again, T is the true kinetic energy of an interacting system and Tg is the kinetic 
energy in the non-interacting system, which we can calculate. All interaction 
effects except the classical coulomb interaction are summed up in the £xc term, 
of which no exact form is known yet. 

Thus, to minimize the energy, the orbitals have to fulfil the following equation 



|18 p.45] 



/ 



]V^^ 



r^dri+yxc(n)-y: 

J ^2 Y 



M 



A 

riA 



\ 



(Pi = Eicpi. (2.20) 



J/ 



In the square brackets, the first term is the electronic coulomb potential, the 
second is the potential that corresponds to the exchange-correlation energy 
Vxc = -^ and the third is the nuclear coulomb potential. These are often 
referred to as an effective potential V^eff^ 



M 



Veff- f^-^drUVxc{r-l)-f^. (2.21) 

J ri2 ^ Tip, 

Up to this point, the Kohn-Sham Theory is exact, there are no approximations, 
except Born-Oppenheimer, so far. If the true analytic forms of all functionals 
were known, the Kohn-Sham scheme would provide the exact energy. The 
approximation is introduced with the explicit form for the unknown functional 
of the exchange-correlation energy. 

2.3 Exchange-Correlation Functionals 

Exchange-correlation functionals are commonly divided into three groups: lo- 
cal (spin) density approximation (LDA), generalized gradients approximation 
(GGA) and hybrid-functionals. The basis for most functionals is the uniform 
electron gas, which is the model system for LDA, non local effects are then added 
in by GGA functionals. 

In the model of an uniform electron gas, electrons move on a uniform positive 
background charge and the electron density is constant. Despite its simplicity 
it is a valuable model because it is the only system where the exact exchange 
functional is known and accurate numerical approximations to the correlation 
functional are available. 

In the local density approximation, Exc depends only on the local electron 
density: 



I P(j^e(p(j^) 



£xc = P(r>(p(r))df (2.22) 
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where £{p{r)) is the exchange-correlation energy per particle of a uniform elec- 
tron gas of density p{r), which is weighted by the probability p{r) of actually 
finding an electron at r. The exchange energy of an uniform electron gas has 



been derived by Bloch and Dirac 1 22 1. It is given by 



3 3/3p(r) 



(2.23) 



The correlation part has been evaluated numerically by Ceperly and Alder 



1 23 1, and analytical expressions have been interpolated. Commonly used LDA 
functionals are VWN, due to Vosko, Wilk and Nusair [24| and the one due to 
Perdew and Wang [25|. Extending LDA to the spin unrestricted case results in 
the local spin density approximation (LSD A), where 



pLSD 



I p{r)£{pa{r),pii{r))dr. 



(2.24) 



Note that in Equations 2.22 and 2.24 the non uniform density p is inserted, but 



e is always from the uniform electron gas. This is a drastic approximation, but 
performs surprisingly well, for example it yields 36 kcal/mol average unsigned 
deviation on the atomization energies of the G2 test-set (compare Hartree-Fock 
78 kcal/mol) |26| . LDA tends to overestimate binding energies. 

In the generalized gradients approximation, LDA is interpreted as the first 
term of a Taylor-expansion, which is continued to the next term as in 



-LDA 



/^ 



Vp 



£-- = £--+ Cxc(p)-^dr. 



P 



2/3 



(2.25) 



which is called gradient expansion correlation. This tends to yield results worse 
than those obtained from LDA, because the density depletion around an electron 
due to exchange and correlation interactions (the exchange-correlation holaj) 
does not fulfil its physical properties. These properties can be enforced by set- 
ting parts of the exchange hole to zero if they are non-negative and by truncating 
the exchange and correlation holes to the correct sum behaviour (one and zero, 
respectively), which is why the these functionals are then said to work within the 
generalized gradient approximation (GGA). This approach performs reasonably 
well ||27| |28[ |29[ |26|. Popular choices of GGA functionals are BP, with Becke's 
gradient correction to exchange |30| and Perdew and Wang's gradient correc- 
tion to correlation |25|, PW91, where both exchange and correlation gradient 
correction are from Perdew and Wang [l25|, and PBE, due to Perdew, Burke and 
Emzerhof |[3l|. 



^hole functions are an essential concept of DFT but would go beyond the scope of this work 
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Additionally, so called hybrid functionals are available, where the total energy 
is mixed from separate results of different levels of theory. The exchange energy 
is mixed from Hartree-Fock, LDA and GGA calculations and the correlation 
energy is calculated from an adopted GGA functional. Mixing parameters are 
fitted to sets of reference molecules. It has been shown, that these function- 
als perform exceptionally well on molecular systems (see references on GGA 
performance above), but they come at the price of usually twice the computer 
time. 

2.4 The Self-Consistent-Field Procedure 

Because the effective potential depends on the density - and thus on the orbitals - 
one has to insert a guess density, solve for the orbitals and re-insert the resulting 
density. This loop has to be continued until self -consistency is reached, i.e. until 
the change in the density matrix or resulting energy between two steps is low 
enough to satisfy one's convergence criteria. This is called the self consistent 
field (SCF) procedure. It has to be stressed here, that convergence criteria in 
practical calculations are only the error bar within the applied theory. Current 
density functionals are expected to yield absolute energies within an error of 
±25 - 50 kcal/mol, relative energies can be assumed to have an error of about 
±5 kcal/mol. Additional errors emerge, when not all calculation parameters 
are fully converged (usually, to save computer time). A prominent example 
is the basis set employed to generate the density. For structural exploration 
a basis set of single functions for core electrons and two functions for valence 
electrons and additional polarization functions (split valence, for example def2- 
SVP from the Ahlrichs set) might be sufficient. But for energetic evaluations at 
least a basis set of three functions per electron plus polarization functions for the 
valence electrons (def2-TZVP) should be used. A quadruple zeta plus valence 
polarization set (QZVP) is considered close to the basis set limit for DFT. 



Several reviews on convergence issues [321 and functional performance |33 



27 28 29 26 1 are available. 



2.5 Geometry Optimization and Transition State 
Searches 

A structure usually has to be optimized after modelling. Optimization means 
to find a stable structure that represents a minimum on the potential energy 
surface. Such optimizations are founded either on numerical derivatives or 
analytically evaluated gradients (and Hessians) of the total energy with respect 
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to the atomic coordinates. Minima are then found with a usual Newton-Rhapson 
or Conjugate Gradients method. 

Finding the transition state of a reaction is a very similar task to a geometry 
optimization, because transition states are saddle points on the potential energy 
surface. One needs to minimize the energy but find a maximum along one 
direction (the direction of the - a priori unknown - reaction path). Practically, a 
Newton-Rhapson like minimum search is performed, with the restriction that 
the energy be maximized with respect to the eigenvector that belongs to the 
(single) negative eigenvalue of the energy Hessian. This is called eigenvector- 
following. 

This technique requires a reasonable guess structure for the transition state 
which has to be close enough to the true TS to exhibit a negative Hessian 
eigenvalue with corresponding (reaction path) eigenvector. A starting point for 
an eigenvector-following search might be a guessed TS structure obtained from 



synchronous transit methods as proposed by Halgren and Lipscomp in |34| 



Here, all atomic coordinates are moved synchronously from their reactant to 



their product position in a linear interpolation scheme. See also Sections 3.1 and 



2.6 Kinetics 

Reaction kinetics are evaluated in terms of transition state theory, which is 
expected to give a good approximation in the case of thermal equilibrium (fast 
energy exchange of a molecule with the surrounding system). A more complex 
theory like Rice-Ramsperger-Kassel-Marcus is not employed because the large 
uncertainties in reaction barriers from DFT calculations would still render the 
results inaccurate. It must be noted here that transition state theory gives an 
upper bound to the rate constantsF] 

For illustration, consider a simple uni-molecular reaction from species A to B 
with one transition state, called X*: 

A ^ X* ^ B 

The equilibrium constant for the first part of the reaction K^ (from the reactant 
structures to the transition state) is 

c{A) Za 



^This section is an excerpt of previous work [3] and only included for completeness 
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where c*, Ca are the concentrations of the transition state (TS) and the reactant, 
z*, Za their partition functions, respectively, and AEq is the difference of ground 
state energies between reactant and transition state. R is the ideal gas constant. 
z^ may now be written as 



z* = z^ 



I — 0-hv/k^T 



where v is the frequency of the vibrational mode of the TS along the reaction 
coordinate and z^ is the partition function without that vibration. Since the 
corresponding force constant is very low, the frequency will also be very low 
and the exponential function is expanded into a series and all but the linear term 
are neglected 



2* - Ztq _,,v/fr„T ^ ^r. n ETTTiT^ - ^t 



hT 



The TS has 2v chances per unit time to dissociate (v in each direction). Due to the 
very low force constant, one may assume, that every chance for dissociation will 
be used. However, only half of the incidents will lead to the reaction products, 
the others will lead back to the reactants. Thus, the rate constant is 



= 2v--K* 
2 

/CB-' 1 -EJRT 

= v -Zf e "' 

hv Za 

h Za 

= ^ . gAS(T)/R . g-AH„ib(T)/RT . ^-EJRT ^2.26) 

h 



which is known as the Eyring Equation. Here, S is the entropy, Hvib the 
enthalpy of the system due to vibrations (including ZPVE) and Eg the activation 
energy of the reaction. AHvib is defined as 

AHvib = Ha{T) - Ea 
and can be calculated by most electronic structure programs in an harmonic 
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oscillator approximation. The energy levels of one vibration are determined by 

en-fico[n + 



2 
2nV [I 

mim2 



nil + ^2 

where co is the vibrational frequency, k is the force constant and [j. is the reduced 
mass. 
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3 Computational Details 

Structure optimizations and energy evaluations were performed with the Tur- 



bomole pSi and DMoP |36 37| programs. Within Turbomole, the b-p functional 



was used. It contains Slater and Dirac's LDA exchange functional [22 38 1 with 



Becke's 1988 gradient correction |30| and Vosko, Wilk and Nusair's LDA cor- 
relation functional [24] with Perdew's 1986 gradient correction [39 1. Molecular 
orbitals were expanded in a basis set of atomic orbitals, in Turbomole a def2-SVP 
(Split Valence plus Polarization) basis was used for structural exploration and a 
def2-TZVP (Triple Zeta Valence plus Polarization) basis was used for quantita- 
tive description. The multipole accelerated resolution of identity approximation 
for coulomb interaction (MARI-/) [40 1 was employed with corresponding aux- 
iliary basis sets ['^Tlj. 

Using DMoP, the BP functional was employed, which has the same exchange 
term as b-p in Turbomole but implements Perdew and Wang's 1992 gradient 



correction [25| to the correlation term. DMoP provides numerical orbitals, of 
which the DNP (Double Numerical plus Polarization) set was used. 

For exploration of the potential energy surface in geometry optimizations and 
transition state searches total energies were converged to 10"*^ auHper SCF cycle 
and to 2 • 10"^ au per geometry step, geometries were relaxed to a maximum 
force of 4 • 10"'^ au and a maximum displacement per iteration of 5 ■ 10""' au. 

In order to evaluate the vibrational modes, the structures were re-optimized 
to converge with a maximum force of 1 • 10""* au and a maximum displacement 
per iteration of 1 • 10"^ au with a refined total energy convergence of 10"^ au and 
10"^ au per SCF and geometry step, respectively. In order to calculate vibrations, 
numeric derivatives of the total energy with respect to atomic coordinates of 
the adsorbate were calculated by means of finite differences. Each atom was 
displaced by 0.05 A. The def2-TZVP basis set was employed. 

A full geometry optimization of the smaller cluster structures with DMoP 
takes 24-48 hours on 16 Opteron 8350 CPUs with shared memory, depending 
on the quality of the initial structure. For Turbomole, the timing is about 72 
hours for both geometry optimization and evaluation of vibrational modes on 
16 Opteron 2218 CPUs with distributed memory. 

^Atomic units (au) are hartree (£h, energy, 4.3597482(26) • 10^^ J), bohr (flo, length, 
5.29177249(24)xl0" m) and hartree per bohr (force), see [42| 
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To find a transition state, three major algorithms are commonly in use: eigen- 
vector following, synchronous transit and nudged-elastic-band. Synchronous 
transit and eigenvector following methods have been applied in this work and 
shall be introduced briefly 

3.1 Synchronous Transit Schemes 

The potential energy surface for a reaction may be coarsely screened by linear 
synchronous transit methods (LST). Reactant and product structures are guessed 
and relaxed. An LST algorithm then finds the maximum energy structure in a 
linear interpolation of the coordinates between reactant and product. Coordi- 
nates are moved synchronously (with the same interpolation parameter). An 
estimation for the Transition State can be found in a series of Conjugate Gradient 
(CG) and synchronous transit steps where the minimum in the CG step is used 
as one image in a subsequent synchronous transit step with quadratic interpo- 
lation between the three points. This is referred to as quadratic synchronous 
transit (QST). Coincidence of the CG minimum and the QST maximum is then a 
convergence criteria for the transition state search, because they provide lower 
and upper bounds to the transitions state energy. 

Obviously, within this theory, one has to supply a priori knowledge or chemical 
intuition to guess possible reactions and model reaction products. In a sense, 
one forces the system to do a reaction, and asks how probable it will be in nature. 

3.2 Transition State Searches using Eigenvector 
Following 

Another technique to find possible reaction paths of a system is to evaluate 
vibrational normal modes. Each mode has an associated displacement vector 
in the orthogonalized coordinate system of the vibrations, its eigenvector. Dis- 
placing the molecular geometry along this vector yields configurations that the 
system will eventually occupy when it vibrates at a finite temperature. One ap- 
proach to find reaction paths is thus to identify soft vibrations that might lead to 
a reaction and calculate the system energy along the corresponding eigenvectors 
{line scanning). When the energy curve along the displacement coordinate has a 
negative curvature, then another vibrational analysis of the distorted structure 
will exhibit at least one imaginary mode. One can then follow this mode to the 



transition state in a constrained geometry optimization as explained in 2.5 

This way, an ah initio scan of the potential energy surface in the vicinity of a 
local minimum can be carried out. The most probable reactions that the system 
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will do are investigated and the energetic evaluation of the transition states is a 
by-product. 

Thus, theoretically, no a priori knowledge or chemical intuition is needed 
to study possible reactions. However, in practice, it takes significant effort and 
human interaction to select meaningful vibrational modes and drive them to the 
actual transition state that leads to a reaction. In addition, this method becomes 
problematic when angular vibrations are expressed in Cartesian coordinates and 
need to be deflected far from the equilibrium structure because atoms will be 
moved on the tangent to the bond angle instead of the real vibrations' coordinate. 
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4 Model System 



Wang et al showed that Copper(II)oxide (CuO) reduces to metallic copper via 
the intermediate formation of Copper(I)oxide (CU2O) when a limited supply 
of carbon monoxide is used as reduction agent [9|. They also mention that the 
reduction of CU2O to Cu is much harder than the reduction of CuO to CU2O. The 
X-ray photoelectron spectroscopy (XPS) data of the ALD-samples from Wachtler 
et al. exhibits a much stronger Cu(I) peak than a Cu(II) one[l|. Therefore, and 
because it is the rate limiting step of the total reduction, CU2O has been chosen 
as a model system and all reactions are modelled on CU2O surfaces. 

Copper(I)oxide has a cubic structure, that was described by Hafner and Nagel 



|43| and Kirfel and Eichhom [44]. The (111) and (100) surfaces exist naturally. 



but the (100) surface is polar and literature is not consistent on the surface 



reconstruction 1 12 , 15 1. 



Thus, for modelling, different Cu20(lll) clusters were used. See figure |4] 
for the Cu20(lll) structure. The Cu20(lll) surface consists of copper-layers, 
which are embedded in two slightly shifted oxygen layers. In this work these 



three layers together are considered one "atom layer", contrary to e.g. |45|. On 
the surface, coordinatively unsaturated copper atoms exist, these are the ones 
with no visible bonds in figure |4|la, they have one coordinative bond to the 
underlying oxygen layer and are surrounded by a hexagon of saturated copper 
atoms, three oxygen atoms in a slightly raised layer and three oxygen atoms in 
a slightly depressed layer. Coordinatively unsaturated copper atoms are called 
Cucus/ oxygen in the raised layer Osuf and oxygen in the lower layer Osub- All 
clusters were stoichiometric and 3 atom layers deep, two with total formulas 
CU70O35 (see figure |4J2 and|4|3), one CU112O56 (not shown) and one CU124O62 (see 
figure |4J4) were modelled. The two smallest clusters differ in the arrangement 
of the outer copper atoms. Common to all clusters are a Cucus species on one 
side and a central Osuf species on the other. These are the expected active sites 
on a clean surface, it is therefore important that they are located centrally to 
minimize border effects. All CU2O surface structures discussed in literature 
can easily be prepared using the existing cluster structures by removing the 
respective Cucus and Osuf species. 
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O(iii) 



b) i (111) 




Figure 4.1: Top (a) and side (b) view of the Cu20(lll) surface, a) shows only one 
atom layer, for clarity; b) shows three atom layers as in the cluster 
structures that were employed. Red atoms are oxygen, orange ones 
are copper. 





Figure 4.2: Top (a) and side (b) view of the first CU70O35 cluster (cluster 1), 

relaxed structure, grey atoms in a) were fixed during the geometry 
relaxation. 
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Cf"*X 



Figure 4.3: Top (a) and side (b) view of the second CU70O35 cluster (cluster 2), 
relaxed structure, grey atoms in a) were fixed during the geometry 
relaxation. 




b) 




Figure 4.4: Top (a) and side (b) view of the CU124O62 cluster, relaxed structure. 
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Another option for a model system would be a periodic slab. In this model a 
unit cell with a few atom layers of the surface is created and periodic boundary 
conditions are applied. A thick vacuum space is left between the surface atoms 
and the top face of the unit cell, see Figure |4J5 for an illustration. The advantage 
is, that no border effects influence the surface atoms because they reside in the 
potential of the infinite repetition of the unit cell. However, sufficiently large 
unit cells have to be used in order to minimize interactions between adsorbed 
species in directions parallel to the surface and interactions between periodic 
surface images in directions perpendicular to the surface. Such systems are 
best described by a plane wave ansatz for the wave function or density in a 
Fourier transformed space. Thus, empty space has to be treated with the same 
amount of basis functions, which is computationally expensive, furthermore 
only few density functional programs are capable of calculating modern hybrid 
functionals from quantum chemistry using a plane wave basis. 

In comparison, the cluster model allows the application of fast quantum 
chemistry methods with atom centred basis sets and the full variety of available 
functionals. Computational effort is usually lower than with periodic systems. 
Energies of reactions on the surface are free of interaction with periodic images of 
the system but are sensitive to border effects introduced with the finite surface 
representation. Thus again, a sufficiently large model has to be created in 
order to minimize such effects. A cluster model corresponds to a low surface 
coverage of the adsorbed species because one adsorbed molecule is subject 
to the investigation while a cluster model represents a high surface coverage, 
depending on the size of the unit cell, an adsorbed molecule may already cover 
all active sites and will always interact with neighbouring cells. 



H 



5?ISi^f' 



ly^ 



Figure 4.5: A slab model for the CU2O surface with adsorbed formic acid, 
vacuum space is cut off for better scaling. 
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5 Results and Discussion 



5.1 Geometry of the CU2O cluster structures 

A geometry optimization will give a first hint, whether a model is reasonable. 
If the structure does not distort heavily and atoms stay close to their bulk (or 
surface reconstructed) positions, it is likely that the model represents a valid 
reproduction of natural surfaces. In a free geometry relaxation of the cluster 
structure, outermost copper and oxygen atoms were displaced by more than 
0.9 A from their bulk positions in directions parallel to the surface. Inner atoms 
moved only slightly, especially the hexagon ring around the central Cucus atom 
was stable. In consequence, to avoid border effects and put inner atoms into 
a more natural electrostatic embedding, shell atoms were kept fix at their bulk 
positions in all subsequent calculations (see figure |4}2a, fixed atoms are greyed 
out). Cucus species moved off their bulk positions slightly when no symmetry 
point group was enforced, which is considered due to a numerical error. 

Only the two smallest of the cluster structures were used throughout this 
work. Relaxed structures of the larger clusters are available, but no reactions on 
their surface were modelled. 

Also several slab models were created with a bulk terminated surface of 
three atom layers. No significant reconstruction took place during geometry 
relaxations of the topmost layer. Again, all terminations discussed in literature 
may easily be created, starting from this model. These models may be used for 
further validation in ongoing studies. 



5.2 Adsorption of formic acid 

Adsorption of formic acid will be the first step towards any reaction mechanism 
on the surface. In order to predict adsorption structures and energies several 
geometry optimization procedures of formic acid on the CU2O (111) surface 
were performed. This procedure yielded two stable adsorption structures. Two 
(A and B) exhibit a large adsorption energy (50 to 85 kcal/mol) at both levels 
of theory and are thus considered to be stable. Both have a bridged structure 
where either the acid hydrogen or the carbon-bound hydrogen binds to one 
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DMoP, BP, 
DNP basis 


Turbomole, b-p, SVP basis 




Cluster 1 


Cluster 1 


Cluster 2 


Structure A 


30 


85 


65 


Structure B 


17 


68 


50 


Structure C 


6 


63 


40 


Structure D 


5 


not stable 


not stable 



Table 5.1: adsorption energies in kcal mol ^. 



OsuF species and the doubly bound oxygen coordinates with a Cucus species. 
Structure C has a high adsorption energy (about 50 kcal/mol) within Turbomole's 
b-p, SVP level of theory but DMoP's BP/DNP level of theory predicts only 6 
kcal/mol adsorption energy. Structure C is not stable within the b-p/SVP level 
of theory in Turbomole, a geometry optimization proceeds towards Structure A. 
The reason for the qualitative discrepancy between the two levels of theory in 



describing structures C and D has yet to be investigated. See Figure 5.1 for the 
adsorption structures on cluster 1, where only parts of the surface are displayed 



for clarity. Table 5.1 lists adsorption energies for the different structures within 
both theories and on the two CU70O35 cluster models. While adsorption energies 
are generally larger within b-p/SVP; both levels of theory agree on the general 
trend that the carboxyl-bridged structure is energetically more favourable. Note, 
that adsorption energy is the negative reaction energy of the adsorption reaction. 
Thus, a positive adsorption energy indicates a stable adsorption. 

The two CU70O35 cluster models have been used to find adsorption structures 
and the results obtained with both agree on structure and the general trend of 
the adsorption energy of the structures. 
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B 



tA^ 



D 



Figure 5.1: Adsorption structures of formic acid on cluster 1, red atoms are 
oxygen, orange ones copper, white ones hydrogen and grey ones 
carbon. 

A acid hydrogen binds to one Osuf species, doubly bound oxygen 
coordinates with a Cucus species. 

B carbon-bound hydrogen binds to one Osuf species, doubly bound 
oxygen coordinates with a Cucus species. 

C carbon-bound hydrogen binds to one Osuf species, acid oxygen 
coordinates with a Cucus species. 
D carbon coordinates with an Osuf species (top view). 
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5.3 Decomposition and Reaction Patlis 



5.3.1 Vibrational Analysis of the adsorbed Formic Acid 
Molecule 

A vibrational analysis of the adsorbed formic acid molecule has been carried 
out. Most importantly, vibrational data is a link to experimental work as it 
allows experimentalists to search for surface species when theoretical spectra 
are available and it allows to verify the theoretical model when spectra of known 
surface species are available. The spectra of the two adsorbed species allow to 
distinguish between them in an experimental spectrum. The normal modes 
from such calculations are also starting points for transition state searches in the 



eigenvector following scheme, as explained in 3.2 
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— freeHCOOH 

— adsorption structure A 

— adsorption structure B 
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3000 



4000 



Figure 5.2: Infrared (vibrational) spectra of formic acid, comparing the free 
molecule with the adsorbed species. 



Figure 5.2 shows the vibrational spectrum of formic acid for comparison of the 
free molecule with the adsorbed structures A and B. For the adsorbed species, 
only vibrations of the adsorbate atoms were calculated, freezing the surface 
cluster's atoms. Mostly, the frequencies change only little, but the trend is to 
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red-shift (soften) for bond-stretching modes and to blue-shift (harden) for H- 
wagging modes. The O-H stretching vibration of the first adsorption structure 
disappears at 3578 cm"^ but a mode at 314 cm"^ appears. It is still a stretch of 
the O-H bond but when the formic acid molecule is adsorbed, the acid proton 
binds to the surface and the formate rest vibrates vertically in this mode. Thus, 
the mass of the vibrating element (49 u) is much larger than that of a single 
proton in the O-H stretch of the free formic acid. The acid proton is likely in a 
soft double well potential of the two neighbouring oxygen atoms, which would 
account for a great part of that effect. 

Additional modes in the range from 40 to 275 cm"^ emerge when the molecule 
adsorbs, but these are translations and rotations of the molecule on the surface 
and do not lead to configurations relevant for reactions. 



5.3.2 Reaction Modelling using Linear Synchronous Transit 

In order to evaluate the model's applicability for the prediction of reaction mech- 
anisms and energies and to gain a first understanding about possible reactions 



on the surface, some of the reactions from section 1.3 were modelled in a syn- 
chronous transit scheme as implemented the DMoP program. It must be noted 
that the results are upper bounds to the reaction barriers but a general trend on 
kinetic probability of each reaction will be discernible. Two of the suspected de- 



composition reactions (1.6 and 1.8) were modelled. Furthermore, in a technical 
application the copper oxide surface will never be clean. Adsorbed hydrogen, 
hydroxyl and water species from wet oxidation pulses as well as organic disso- 
ciation products of the copper precursor are likely to saturate the surface. Thus, 
some of the reactions with previously adsorbed species have been modelled 



too. These are reactions |1.12t p34l and prT6l Note that this is not intended as an 



exhaustive search for a reaction mechanism, but as a starting point for ongoing 
research and to gain an initial understanding on how reactions may be modelled 



with the surface structure at hand. Table 5.2 lists reaction energies and barriers 



as predicted by an LST-I-QST search and estimated rate constants k (see below). 



All modelled reactions are illustrated in figure 5.3 

Most of the reactions are exothermic, but the pure acid deprotonation to yield 
an adsorbed proton is slightly endothermic. In contrast to the gas phase reaction, 
the decomposition to CO and H2O is exothermic but still has a high activation 
barrier of ±80kcal mol" . Both the acid deprotonation to form molecular hy- 
drogen with an adsorbed proton and to form water with surface oxygen and an 
adsorbed proton are strongly exothermic but hindered by high reaction barri- 



ers. Reaction |1.16| is an actual reduction reaction, in the sense that it removes an 



oxygen atom from the surface. A reaction with adsorbed OH to form water is 



31 



slightly exothermic, the reaction barrier has yet to be calculated for this reaction. 
An accurate prediction of reaction kinetics from first principles is difficult. It 
is not possible to report actual error bars on any electronic structure calcula- 
tion, and even when assuming a very small error in the reaction barrier of 1 
kcal/mol the reaction rate constant has an uncertainty of one order of magnitude 
because the activation energy enters an exponential function. However, from 
the available data, a rough estimate of reaction rates can be made according to 



the EYRiNG-Equation (2.26 ), assuming zero for AS and AHvib- These are expected 



to introduce less than ±5 kcal mol in the exponent from the findings in |p 



It can be predicted, that the pure deprotonation by reaction 1.6 will not be 
limited from the intrinsic reaction rate but rather from the amount of free ac- 
tive sites on the surface. In contrast, the latter two deprotonation reactions are 
certainly kinetically unfavourable because of their high activation barrier. Re- 



action 1.16 might proceed slowly, a rough estimate for the reaction rate at 400 K 
is 10-1^ 



s\ But reaction 



1.12 



has a far too high barrier and will not proceed at 
400 K before other reactions take place. A reaction rate constant is estimated to 
about 10-33 g-i 

To present a more comprehensive quantity, the theoretical half -life of a surface 
population with one of the reactants that undergoes only one reaction is given by 
f 1 = In 2/fc. The half-life times for the three high-barrier reactions are lO^*^ s, lO^^ s 

2 ° 

and 10^^ s, so even if the reaction barriers are overestimated by 10 kcal mol" , no 
significant conversion will take place through these reactions. 

The results show, that a kinetic evaluation of model reactions within the 
limitations of the applied theory is possible with the model system at hand and 
can be used to test for probable and improbable reactions and mechanisms. 

Concerning a reduction mechanism for the copper oxide layer, it is clear that 
other mechanisms with lower kinetic hindering must be existent or catalytic 
influence will be necessary in a technical application. 



# 


reaction 


AE 


Ea 


k 


1.6 


HCOOH(,,,) -^ H^^,^^ + HCOO(,,3) 


7 


10 


W 


1.8 


HCOOH(,,3) -^ CO(,,3) + HP(^,^^ 


-16 


80 


io-3i 


1.12 


HCOOH(,,, + H(^,^) -^ HCOO(,,3) + H.^^^ 


-35 


84 


10-3^ 


1.14 


HCOOH(,,3) + OH(,,3) -^ HCOO(,,,) + Hp^^.^^ 


-3 


n.a. 


n.a. 


1.16 


Og^p-HOOCH + H(^,^) -^ HCOO(,,3) + Up^^^ 


-18 


47 


10-14 



Table 5.2: reaction energies and barriers in kcal mol ^ and rate constants at 400K 
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Figure 5.3: Overview of the modelled reactions, with reaction energies and bar- 
riers in kcal mol"^ and reaction rate constants in s"^ 
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All reactions were modelled using uncharged structures. Mulliken partial 



charges of some of the adsorbed species are listed here. After reaction 1.6 the 



adsorbed hydrogen has a charge of 0.286 e and the formate rest has a total charge 
of -0.472 e, where e is the elementary charge. The adsorbed single hydrogen in 



the reactant structure for reactions 1.12 and 1.16 has a charge of 0.171 e. In the 



reactant structure of reaction 1.14 the hydroxyl oxygen has a charge of -0.566 e 
and the hydrogen has a formal charge of 0.264 e. A thorough discussion of 
the bonding mechanisms and electron density distributions should be done in 
future projects but would go beyond the scope of this work. 



5.3.3 Transition State Searches using Eigenvector Following 

Line scans of the vibrational modes of adsorption structure 1 have been per- 
formed. 

A complete assessment of the chemically interesting modes was not possible 
in the time frame of this work. One transition state between two adsorption sites 
could be localized. However, no converged structure of transition states for de- 
composition reactions can be reported from eigenvector following calculations 
of the negatively curved line scans. 

The results show, that a full evaluation of reaction barriers from an eigenvector 
following scheme will take significant effort that would go beyond the scope 
of this work. Eigenvector following methods might however prove useful for 
further validation and refinement of the results obtained from synchronous 
transit methods. 



Figure 5.4 shows the transition state for mobility of the adsorbate from one 



Cucus to another. The oxygen atoms in the formic acid will change roles during 
this reaction. The reaction barrier is 24 kcal mol"^. 




Figure 5.4: Transition state for mobility of the adsorbate from one Cucus to 
another. 
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6 Summary and Outlook 



Four cluster models for a copper(I)oxide (111) surface have been designed, of 
which three were studied with respect to their applicability in density func- 
tional calculations in the general gradient approximation. The two smallest 
systems have been found to be computationally feasible and produce quali- 
tatively matching results, investigation of the other structures is an ongoing 
project. Further validation will be necessary for quantitative results. The larger 
clusters are similar in design to the smaller ones and may by comparison pro- 
vide information about border effects and convergence with system size. The 
models exhibit the active sites of the Cu20(lll) surface and are capable of pre- 
dictions about adsorption and subsequent reactions. Reaction energies and 
barriers can be estimated in the scope of density functional theory calculations. 
Further validation against a slab model with periodic boundary conditions and 
experimental data should be performed in further work. 

Formic acid adsorption on these systems was modelled and yielded four 
different adsorption structures, of which two were found to have a high (50-85 
kcal/mol) adsorption energy. Vibrational spectra of adsorption structures A and 
B are available. These are important for experimental distinction of the two 
structures and allow further validation of the calculations once experimental 
data is available. 

The energetically most favourable adsorption structure (A) was further inves- 
tigated with respect to its decomposition and a few reactions with adsorbed H 
and OH species using synchronous transit methods to estimate reaction barri- 
ers and single point energy calculations for the reaction energy. First, the acid 
deprotonation of formic acid on the surface was modelled. The reaction is re- 
ported to be endothermic (A£ = 7 kcal mol" ) with a modest reaction barrier 
of 10 kcal mol" . It is thus not limited by reaction kinetics but by formic acid 
concentration and adsorption sites. The decomposition to CO and H^O was 
found to be exothermic (AE = -16 kcal mol" ) but is inhibited by a high barrier 
of 84 kcal mol" . The reaction with an adsorbed proton to form molecular hy- 
drogen causes an exothermic energy change of -35 kcal mol" but it also has a 
high barrier of 84 kcal mol" . One actual reduction reaction that would remove 
a surface oxygen to form water with the acid hydrogen and another adsorbed 
proton is also exothermic (AE = -18 kcal mol" ) but hindered by a barrier of 
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47kcalmol \ 

Although reaction barriers, as calculated so far, are upper bounds and rate 
constants are strongly dependant on the barrier, no change in the qualitative 
trend is expected because the barriers are very high. Thus, a reaction mecha- 
nism for the reduction of copper oxide by formic acid must go through different 
reactions not considered here or needs catalytic support in the surface, for ex- 
ample by ruthenium inclusion as investigated in fSf. 

The results show, that the presented cluster models can be employed for 
electronic structure calculations at the density functional level to compute reac- 
tion energies and barriers. All Cu20(lll) surface reconstructions discussed in 
literature can easily be modelled using the existing structures. 

Thus, the present work may serve as a basis for detailed studies of reactions on 
the CU2O (111) surface which might include further reactions involving different 
surface species, more spectroscopic data for comparison with experimental work 
and a thorough investigation of the chemical bonding mechanisms. If necessary 
more accurate electronic structure methods than DFT might be employed. 
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